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Introduction 
Icing Observations 


r 


• Many engine power-loss events reported since the 1 990’s 

• Mason et al. 1 hypothesized how power-loss events can result from ice 
crystals entering the engine core 

• Ingestion of ice into engine is studied 

• NASA conducted tests at NRC’s RATFac 

• Observed environmental conditions changed with cloud activation 

- Air temperature change 

- Air humidity change 

- Particle water content change 

• Hypothesis : Thermal interaction between air and cloud particles 


1 Mason, J., Strap, J., and Chow, P., “The Ice particle Threat to Engines in Flight,” presented at 44 th AIAA Aerospace Sciences 
Meeting and Exhibit 2006, USA, January 9-12, doi: 10.2514/6.2006-206 
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Thermal Model General Description 


• Model couples air and cloud particle conservation equations 

- Mass, energy and momentum 

• Simulates icing tunnel (applicable to engines too) 

• Ice, water, and humid air mass broken into fundamental CV units 

- Uniform particle size (MVD representation) 

- Full particle size distribution 



Single particle depiction (uniform particle size cloud) 



Full particle size distribution particle depiction (multiple particle size cloud) 
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Assumptions 


• Air and water vapor are treated as ideal gases 

• Air is continually well mixed 

• No supersaturation 

• 1 -D air and particle flow 

• Dilute system (no particle interaction) 

• Particles are spherical 

• Discrete particle size distribution (bins) 

• Uniform temperature within the particle 

• Supercooling can occur 

• Mixed phase particles are not spatially resolved 

• Phase change occurs at particle surface at particle temperature 

• Adiabatic tunnel walls 

• The flow of particles and air is a continuous stream 

• The fundamental CV is adiabatic and mass is conserved 

- Provisions for heat and mass transfer added using source terms 
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Particle Conservation of Energy 
Single Particle Formulation (1/2) 




Change in the particle’s enthalpy, is due to the convective heat transfer 
and latent energy exchange due to mass transfer 

dH p _ 

o , — Rconv R latent 


• Rate of mass change with a water surface is proportional to the difference 
in vapor pressure between the particle surface and the ambient air 

Ice. PpCp—~ ud h(T air Tp ) + ud h m p a i r L su pi(^co a i r co p ) 


Mix. PpL m eit 


nd 3 dri p 
6 dt 


Ted h(r air Tp) + Ted h m ,P a i r Lsubl/ev a p(.^ a i r ^p) 


Water. PpC p ~ ^7 ?cd h(T a i r Tp) + 7 id h m p a i r L eva p(^co a i r co p ) 
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Particle Conservation of Energy 
Single Particle Formulation (2/2) 


• Empirical heat and mass transfer expressions for flow over a sphere 2 


- JVu = = 2 + 0.6 Re 1/2 Pr 1/3 

P- air 


- Sh = 

- Re = 

- Pr = 

- Sc = 


^=2 + 0.6 Re 1/2 Sc 1/3 

Dab 

P air\ v air ~ v p \d 
Pair 

C air Pair 
kair 
Pair 
Pair Dab 


c Dair ~ 

( Dsurf ~ 


MW : ; 


water 


MW air 
MW 


water 


MW, 


air 


wv,air 
r air 

Dwv,surf 
D air 


2 lncropera, F., and DeWitt, D., “Fundamentals of Heat and Mass Transfer, Fourth Edition,” (New York, John Wiley & Sons, 1996) ISBN: 0471304603. 
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Air Conservation of Energy 
Single Particle Formulation 


• Change in the air enthalpy is due to the convective heat transfer and 
the sensible energy change of the water vapor mass that has 
changed phase 

SHair _ , 

~ Qconv ' Rwv,sens 

• tfwv.sens is thermal mixing 

- "lairCair = nd 2 h(T p - T air ) + C wv OT 

• m wv can be mass from evaporation or mass to condense 
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Conservation of Mass 
Single Particle Formulation 


• Change in particle mass due to vapor phase change 

“ IfT = PpYt Gt~) = nd 2 h mPair(<»air ~ <Wp) 

• Change in air mass is opposite the change in particle mass 

_ = _ 9mp = 7 Td 2 h m p air ((i) p - 0 ) air ) 
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Conservation of Energy 
Particle Distribution Formulation 


• M VD vs particle distribution 

- Smaller particles have greater SA/Vol ratio -> faster transfer response 

- Cumulative differences will add up to a different final result 

- Fraction of the total water mass calculated for every particle size bin 


An energy balance equation for every particle size / 

— Pp,i^p,i~2 f -!l ~ nd-i hi(T air — T p i) 4- 7r di h m ip a i r Li((i) a i r — 


One air energy equation contains the sum of all the particle heat 
transfers and vapor sensible energy transfers. 

' (#;) 


r dT air 
m a ir C air 


dm , 


— Y™ ird?h ( T ■ — T ■ ^ 4- uniwv ’ i C 1 ?- 1 r . f)T 

— Zji=i y LLL i u i\ 1 p,i l airJ ~ Q t JT r :„ U1 


air 


n = number of particle size bins 
#i = number of particles in the i th bin 
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Conservation of Mass 
Particle Distribution Formulation 


• A mass balance equation for every particle size /' 

— —gf— — ltdi h m ip a l r ((l) a i r — 

• One air mass equation contains the sum of all the particle mass 
transfers. 

_ d m air yn 1 \ 

dt ~ ^ =1 L dt 
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Conservation of Momentum 


Conservation of momentum equation is solved in reference to particle 


F Fclrag F Fg — TtlpCL — Pp 
~ Fg = 0 

— F dr a cj ~ 2 Pair^ ACq 
~ U V a ir 


nci 3 dv p 
6 dt 


- A = nd 


24 

- c D =- + 


Re 


■4Si 

{Re A ' 3 

V5.oJ 


+ 


0.41 if — — — ) 

V263,000 / 


-7.94 


1 + 


(-5S-) 

\263,000 / 


- 8.00 


+ 


Re 


0.80 


461,000 


(Morrison correlation 3 ) 


dVp 3 P airC d(v air ^p) 


V 
dt 


p p d 


3 Morrison, F., “An Introduction to Fluid Mechanics,” (New York, Cambridge University Press, 2013) pg. 625, ISBN: 1107003539. 
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Supercooled Freezing Formulation 


r 


• Four stages in supercooled freezing 

1. Sensible liquid cooling (below 273.1 5K) 

2. Latent heat release 

3. Latent freezing (at 273.1 5K) 

4. Sensible ice cooling 



• Homogeneous crystallization temperature as a function of diameter 4 

- T hc = 7.2015 ln(d) + 214.64 


Expressions to determine melt fraction after latent heat release 

— H sens, super iy* Cp(T")dT 


Vp = 1 - 


TTLpLmelt H sens, super 
TtlpLmelt 


4 Heverly, J., “Supercooling and Crystallization,” Transactions of American Geophysical Union 30(2): 205-10, 1949. 
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Algorithm 


r 


• W ritten in MATLAB version R201 4a 

• Solves conservation equations using built-in ODE45 solver 

• Relative and absolute convergence tolerance of 1 0 -12 

• Mass transferred between the air and particle(s) balanced to 1 O' 15 

• Energy transferred between the air and particle(s) balanced to 1 0 -4 

- Model accuracy dependent on accuracy of property values (C p , L heat , etc.) 
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Parametric Analysis 
Baseline 1 Tests 



Test # 

Parameter Chanqed 

Units 

Value 

Baseline 1 

1 

Slip Velocity 

m/s 

25 

2 

Particle Temperature 

K 

273.15 

3 

Pressure 

Pa 

28,000 



o 0 5 10 15 20 

Timp 


Takeaway: Wet-bulb temperature determines state of particle 
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Parametric Analysis 
Baseline 2 Tests 



Units 

Baseline 2 

Particle Diameter 

pm 

10 

IWC 

g/m 3 

1 

Pressure 

Pa 

88,000 

Relative Humidity 

% 50 

Air Temperature 

K 

280.15 

Particle Temperature 

K 271.15 

Slip Velocity 

m/s 

5 


Parameter Changed I Units 


Baseline 2 

4 

IWC 

g/m 3 

5 

5 

Relative Humidity 

% 

80 

6 

Air Temperature 

K 

271.15 


Value 


Test # 



Takeaway: Air thermal mass » Particle thermal mass 
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Thermal Model - RATFac Modifications 
Air Heat, Air Velocity, Air Vapor Mass Sources 


• Air velocity, air temperature, and vapor content vary from injection to 
test section at RATFac (no icing cloud) 

• CFD model previously written to approximate v air and T air 

• Centerline values normalized 

• Water vapor profile analogous to normalized air temperature profile 

• Approximate changing values as sources in expressions (no cloud) 

- Valid for tests Mach ~ 0.2 and P total ~ 44000 Pa at test section 




Normalized air velocity at centerline (without an icing cloud) Normalized air temperature at centerline (without an icing cloud) 
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Model - Experiment Comparison 


|| Units 

Scan 877 

Scan 983 

Scan 1003 

Tair,inj 

K 

262.0 

256.1 

256.9 

Tair, target 

K 

288.4 

278.0 

277.9 

SH inj 

(RHinj) 

Qvapor /kQdry air 
(%) 

0.07 

(2) 

0.07 

(5) 

0.07 

(4) 

SHt ar get 

(RHtarget) 

Qvapor /kQdry air 
(%) 

4.07 

(16.1) 

2.88 

(35.4) 

3.81 

(46.5) 

GWCjpj/target 

g/m 3 

2.1 

2.4 

3.0 

LWCjpj/target 

g/m 3 

1.0 

0 

1.9 

IWCj n j/t ar g e t 

g/m 3 

0 

8.4 

8.6 

Hinj/target 

- 

1.0 

0.0 

0.18 

P 

Pa 

42806 

66478 

65934 

v air,inj/target 

m/s 

86.8 

85.7 

84.1 

Twater,inj 

K 

278.15 

- 

278.15 

T ice,inj 

K 

- 

256.15 

256.15 

MVD wa ter,inj 

gm 

40.0 

- 

40.0 

MVDicejpj 

gm 

- 

45.5 

45.5 



Units 

Scan 877 

Scan 983 

Scan 1003 

AT air , exp 

K 

-2.6 

-2.9 

-4.4 

A"T a j rj S j m 

K 

-0.54 

-0.75 

-0.88 

AGWCexp 

g/m 3 

0.5 

1.0 

1.2 

AGWC sim 

g/m 3 

0.13 

0.16 

0.21 

ALWCexp 

g/m 3 

-0.5 

0.5 

0.2 

ALWC sim 

g/m 3 

-0.13 

0.0 

-0.09 

AlWCexp 

g/m 3 

0.0 

-1.5 

-1.4 

AIWC sim 

g/m 3 

0.0 

-0.16 

-.12 

Anexp 

- 

0.00 

0.07 

0.05 

Arisim 

- 

0.00 

0.00 

-0.004 


• Conditions 

- Scan 887: Water sprayed into low P, low RH, very warm air 

- Scan 983: Ice sprayed into medium P, medium RH, warm air 

- Scan 1003: Water & Ice sprayed into medium P, medium RH, warm air 

Takeaway: Model accounts for -20% of experimentally observed changes 
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Sources for Model - Experiment Discrepancy 


• Experimentally observed water/ice film on tunnel walls 

- Increased residence time for evaporation 

• Poor approximation for the vapor mass source term 

• Non-spherical particles increase heat and mass transfer 

• Turbulence at spray nozzle may enhance transfer 

• Uncertainties with experimentally measured values 

- Independent corroborating measurements minimize this possibility 
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Conclusion 


Thermal model that couples particle and air conservation equations 

• Simulates: 

- All phase change types 

- Supercooled and normal freezing 

- Single particle and full particle distribution sprays (ice, water, combined) 

- Complicated icing tunnels with energy and mass sources 

• Air temperature, pressure (air mass), and RH dominate -> T wb 

• Model compared to experiments conducted at NRC 

- Simulated ~ 20% of the cloud and air changes observed experimentally 

- Reasons for discrepancy are offered 

• Future work to determine sources for discrepancy 

• Model can be modified to simulate other icing facilities (PSL) 
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